Dynamical magnetic excitations of nanostructures from first-principles 
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Abstract 



functions, we devise a real space method to investigate spin dynamics. Our scheme enables one to 
deduce the Coulomb potential which assures a proper Goldstone mode is present. We illustrate with 
application to 3d adatoms and dimers on Cu(lOO). 
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Very recently new experimental methods for the study of magnetic excitations in the 
nanoscale and subnanoscale length regime have appeared. Examples are spin polarized elec- 
tron energy loss spectroscopy (SPEELS) [1] and inelastic scanning tunneling microscopy 
(STM) p, yl]. In both methods, the quantity probed is the transverse susceptibility x that, 
in linear response, describes the amplitude of the transverse spin motion produced by an 
external magnetic field Bext of frequency u. To calculate x is a major computational chal- 
lenge, if one wants to employ time dependent density functional theory (TD-DFT) Q] or 
many-body perturbation theory (MBPT) based on DFtQ]. Thus one sees very few DFT 
based calculations even for bulk systems. Is-Q]. There are studies that use empirical tight 



binding theory (ETB) 



10l-ll3l| based on MBPT. In all schemes, a similar master equation 



must be solved. Its solution maybe be written in schematic notation. 



X = Xo(l - ^Xo)" 



(1) 



In ETB MBPT, Xo is the single particle susceptibility for a mean field ground state, and x 
is the susceptibility generated by the Random Phase Approximation (RPA). In TD-DFT |4|], 
Xo is the Kohn-Sham susceptibility and x. the total susceptibility, is exact in principle if the 
full exchange and correlation kernel U is known. In practice, one often invokes the adiabatic 
local spin density approximation (ALDA) or rarely one uses the RPA approach[9]; U can 
also be viewved as a parameter whose value is in the range of leV/iJ,B for 3d transition 
elements[15]. In the ETB MBPT, U is the effective Coulomb interaction that is the local 
exchange splitting divided by the magnetic moment. Within the ETB where x is generated 
through the RPA the Goldstone theorem is satisfied exactly (the scheme is a conserving 
approximation in many body theory wherein the Ward identity is obeyed) so, for instance, 
zero wave vector spin waves have precisely zero frequency with spin orbit coupling set aside. 
In DFT based methods, in practice the Ward identity is not strictly obeyed so the Goldstone 
theorem is violated in numerical studiess. One then adjusts U in an ad hoc manner so x 
is compatible with the Goldstone theorem. For instance, Sasioglu et al. correct U by 
45% in their study of bulk Ni while Buczek et al. do not calculate it but mention a finite 
frequency shift for the Goldstone mode in the range of 5 — 10 meV js]. This procedure can 



be dangerous 



14l |. and below we introduce a scheme that constrains U so no adjustment are 



needed. 

We report here a new computationally attractive method that allows us to address mag- 
netic excitations without the need to obtain empirical tight binding parameters from other 
electronic structure calculations. We begin with the the Korringa-Kohn-Rostoker single par- 
ticle Green function (KKR-GF) [l6| which contains an ab-initio description of the electronic 
strucutre, and we introduce an Ansatz that allows us to generate a one electron Green func- 
tion whose structure is similar to that which enters the ETB methods. We are then led 
to a computationally tractable and flexible real space scheme for describing spin excitations 
in nanoscale structures. Along the way we derive a criterion for generating an effective U 
compatible with the Goldstone theorem. Other approximate schemes have been introduced 
to address spin excitations within TD-DFT (sl or MBPT based on DFT j^]. Within these 
approaches the calculation of dynamic susceptibilities is still very difficult, and limited to very 
small systems. Our scheme is readily applied to bulk materials, to surfaces with adsorbed 
films, and it is a real space formalism ideal for diverse small nanostrucures. As an initial ap- 
plication, below we explore the spin dynamics of single adatoms and dimers on the Cu(lOO) 
surface. For instance, experimentally, Hirjibehedin et al. [3] employed the semi-insulating 
CuN substrate that differs qualitatively from the more commonly encountered metal sur- 
faces. Their data can be interpreted with a model Hamiltonian which describes an atomic 
like localized moment with integer spin, weakly perturbed by the substrate, whereas strong 
coupling to substrate electronic degrees of freedom qualitatively change the character of the 
spin excitations in systems we explore here: the magnetic moments are neither integer nor 
half integer; our approach addresses systems where the model Hamiltonian used in ref. 
fails qualitatively. 

Our method is based on two assumptions: (i) if the external applied magnetic field is 
constant within each atomic cell of it is spherically symmetric in each cell, one can develop 
a description that maps the calculation onto a scheme discussed many years ago by Lowde 
and Windsor jl?!, though now we have incorporated an ab initio description of the electronic 



12], describes the underlying 



structure into Xo- This structure, as discussed more recently 
physics correctly. With our method the size of matrices involved in solving the Dyson equation 
is modest, thus opening the way to the investigation of large atomic nanosructures with a full 



description of their electronic structure included. Then (ii) In the KKR-GF method Q, the 
Green function is written in terms of products of energy dependent solutions of Schrodinger 
equation in the unit cells; both the regular and irregular solutions enter. We propose here an 
expansion in terms of energy independent d like wave functions we chose to be the regular 
solutions of KKR-GF theory evaluated at the Fermi energy. Our focus is on low energy 
excitations of 3d moments, so this choice is appropriate. 

The first step in solving the Dyson equation is to calculate Xo that is a convolution of 
two Green functions 

X'^{r,r';uj) = j dzf{z)[Gfj{r,r;z + uj)\mGj,{r,r:z) 

+ lmGj(r,r';z)G/(r',r;z-w)] (2) 

where f (z) is the Fermi distribution function, G and G' represent the retarded and advanced 
one particle Green functions connecting atomic sites / and J and ImG = — ^(G — G~). Xo can 
be separated into a sum of two terms: /i which involves Green functions that are analytical 
in the same half complex plane and a non-analytical I2 jl2]- For positive frequencies: 

i^r.r^u) = ^ p dz[f{z-u)Gl{r,r-z)Gl{r,r;z-uj) 

- f{z)G\;{r, r; z)Gj*(r, r^; z - u)]) (3) 

and 

i!i{r, r: uj) = ^ r dzGjjir, r: z + u){Gl{r, r, z) - G];{r, r- z)) (4) 



In order to improve numerical stability, the two terms differ a bit from those presented in 
Refill and they require one less Green function. Such a separation is ^tractive since /i 
can be calculated through use of a regular contour in the complex plane [18] with a modest 
k- and energy-mesh. I2 can be calculated along a line parallel to the real axis that does not 
require much computational efforts since the integration is limited to a small energy regime 
defined by u. 

The green functions are provided by the KKR-GF method |l^: Gjj{f,f';z) — 
-iVzR'^{r<]z)H'^{r>\ z)5,j,Lu + R'^ir; z)GB'^^'{z)RJ^'{f; z) where Ge is the struc- 
tural Green function. Here the regular R and irregular H solutions of the Schrodinger equa- 
tion are energy dependent, and this makes the calculation of Xo in Eq{T]tedious and lengthy. 



Our Ansatz expresses the Green functions in terms of energy independent wave functions 
(p such that: Gjj{f,f';z) ^ J2lli^''~(^^u'~^ where Gy''^(z) is generated from 
r ^Ll r\ r / n^,L;♦/->^ ■ Since the terms in the denominator are normalization factors 

J dr(p"-*{r )(!>"- (f) j dr'(p''-i(r')(p''-i*(r) 

instead of working with (p''-{f) we introduce ip''-{f) = (p''~{f)/{J d f (!>''-*{ f)(p''-{f)) ^ where we 
choose 0{f) = R'^{r, Ep), i.e., the c/-regular solution of the Schrodinger equation. This is 
appropriate for the calculation of the d-block of the susceptibility. 

n 

Within TD-DFT [4], Eq.[T]has been derived via a variational linear response approach. 
Here, a magnetic system, say magnetized initially along the z-direction with a magneti- 
zation mz{f) and a charge density Do(r). is subjected to a small time-dependent external 
transverse magnetic field Bext{r\ f) which induces a transverse magnetization mx,y{r, t). 
The latter quantities are connected to each other through the magnetic response func- 
tion: x'K^^'^^') = 5^xy[Sext](^f)/<56ixt(/^f')leext=o,no- Within the atomic sphere ap- 
proximation (ASA) and assuming a spherical external magnetic field this simplifies to 
X'-'irt, r'f) = 47r(5/77;^[eext](/'t)/(5eixt(r't')|ee,t=o,no. where x'-' is a sum over all angular 
momenta components XIlli X''"'"^'"''"^'"- The same procedure is repeated for the magnetic 
response function Xo of the Kohn-Sham noninteracting system and after a Fourier transform 
with respect to time we obtain a form that resembles the approach of Lowde and Windsor 
which is very often used in the tight-binding simulations of magnetic excitations 
now have 

X'{r,r';u) = x'^{r,r';u) 

+ 5Z / / dr"'Xo{r, r"; u)U"{r", r'"; u)x'\r"' , n u) (5) 

that involves only site dependent matrices. U'\r,r';u) is a functional derivative given by 



m 



Jj. We 



l^M\r^l) \Be.t=o.na that simplifies within ALDA Q to ^^^^^^5rj'5ij assuming an initial 
collinear magnetic configuration; Beff is the magnetic part of the effective Kohn-Sham po- 
tential {V^ff — V^ff). U can be considered as an exchange splitting divided by the magneti- 
zation. 

It is interesting to note that there is another way to determine the correct form of U. 
Indeed, to ensure the realization of the Goldstone mode, it is necessary to have the following 
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sum rule fulfilled at u = |l4| : 



E / dr'x'^i^' ^' 0)Biff{r': 0) = m',{r: 0) (6) 

this expression has been obtained after multiplying both sides of Eq. E] by B^^ff{f;u = 0), 
integrating over f , summing up over all sites j and using the following Dyson equation: 
= Gj, + Gi^BeffG^. Eq. [6] can be rewritten in a matrix notation: 

U' = r-^m', (7) 

with r = Xo/^z at w = and U'{r) = -^^f^rfj^- We then recover the result extracted from 
the ALDA: Eq. [7] provides a means of calculating a value of U which leads to satisfaction of 
the Goldstone theorem. The right hand side requires knowledge of the magnetization and 
the static susceptibility. Stated otherwise, the correct U is the one with the lowest eigenvalue 
of the denominator of Eq{T] associated with the magnetic moments as components of the 
eigenvectors. Such a derivation is advantageous since it would work for more complex systems 
involving non-equivalent atoms. 

Assuming the expansion in terms of energy independent wave functions described previ- 
ously, the final Dyson equation simplifies after some straightforward algebra into a strictly 
site dependent equation 

f = f + f o^f (8) 

where the d-block of the dynamical susceptibility Xo(''''''i^) given by 

^p'^^{r)'ip'^^*{r)xo{oj)'ip^^^*{r')'ipi.^{r') and U can be calculated once for every atom from 
the previous sum rule, Eq. [71 It can be understood as a Stoner parameter and gives once 
more a justification for the approach used by Lowde and Windsor[l7(]: i.e. the effective 



intra-atomic Coulomb interaction is expressed by only one parameter. 

To illustrate our scheme, we shall investigate the magnetic excitations of 3d adatoms and 
dimers deposited on Cu(OOl) surface with the theoretical LDA lattice parameter (6.64 a.u.). 
First we examine the spin dynamics of the magnetic moment bearing adsorbates Cr, Mn, Fe 
and Co. The atoms are positioned in the fourfold hollow sites with the first nearest neighbors 
around every adatom included when evaluating the full real space Green function. 



TABLE I: Comparison between magnetic moments (in /ig) of adatoms calculated from the KKR- 
GF method or following the projection scheme discussed in the text. Values of L/'s (eV/z^e) ^-re 
also shown. 

Cr Mn Fe Co 

Mrf: KKR-GF/model 4.0/4.0 4.1/3.9 3.00/2.7 1.7/1.6 

Mtotai 4.0 4.2 3.0 1.7 

-U 0.90 0.97 0.97 0.97 



As shown in Table [D the ^-contribution to the total moment is the most relevant and 
is nicely reproduced by the projection of the Green functions into our choice of wave func- 
tions. The U's calculated are very close to leV//XB. a value that is v ery often used in the 



literature 



121 ] which was extracted empirically from data by Himpsel 



151]. Himpsel related 



the energy splitting between majority and minority spin states observed in photoemission to 
the magnetic moment of 3d atoms. In Fig. [IJa), we plot the energy splitting between our 
calculated majority and minority spin states as a function of magnetic moment, to find a lin- 
ear relation with slope very close to that in the corresponding figure in Himpsel's paper 
This provides us with confidence in the scheme set forth in this paper. 

We note that many years ago the spin dynamics of a moment embedded in a paramagnetic 



id- 



host was described 



201 ] within a framework conceptually similar to that used here. The local 



response displays a g shifted Zeeman resonance, broadened very substantially by decay of 
the coherent spin precession to particle hole pairs, whereas the total moment of the system 
precesses with g=2 and zero linewidth. Our calculations explore the local response of the 
moment. Local probes such as STM thus examine aspects of spin dynamics inaccessible to 
spectroscopy based on excitation by long wavelength radiation such as microwaves or light; 
these sense the precession of the total moment. 

In Fig. [Hb), we show our calculations of the resonant response of the local moments 
for the four adatoms we have investigated. The resonant frequency scales linearly with the 
applied DC field, as does the width of the structure. For this reason, as in an earlier study [I4], 



we can apply a DC field of sufficient strength to move the resonance up to frequencies that 
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allow numerical study. A simple scaling reduces the results to the regime of physical applied 
fields. The g shift can be positive (Co) or negative (Cr, Mn and Fe). (For the field we have 
applied, a g value of 2 would provide a resonance at 13.6 meV.) These results agree well 
with earlier tight-binding simulations performed on Mn and Co adatom [12]. The width of 
the resonances is controlled by the local density of states [20], and is thus influenced by the 
position of the d levels relative to the Fermi energy. Thus the Co and Fe resonances are 
quite broad, since their minority spin levels intersect the Fermi level, whereas those for Mn 
and Cr are much sharper since for these adatoms the Fermi level lies between the majority 
and minority states. 

We have also explored nearest neighbor dimers of the 3d adatoms discussed above, assum- 
ing a ferromagnetic ground state. For the dimer, there are two resonances, a zero frequency 
acoustical mode, and a high frequency optical mode that is damped by decay to Stoner 
excitations. The position of the optical mode provides information on the stability of the 
assumed ground state. We find the optical mode at negative frequency for Cr and Co. This 
informs us the ferromagnetic ground state is unstable, whereas for the Fe and Mn dimer 
the mode resides at positive frequency so ferromagnetism is stable. This test of ground 
state stability will be useful for studies of more complex structures. We have calculated 
adiabatic exchange integrals J in an effective Hamiltionian we write as —Jei ■ 62, with ei and 



62 unit vectors. 2111 ^Ve find J negative for Cr (J = —19.8 meV) and Co (J = —14.9 meV), 
consistent with instability of ferromagnetism in the ground state, while positive values were 
obtained for Fe (J = 30.4 meV) and Mn (J = 16.3 meV) so we confirm ferromagnetism is 
stable. For positive frequencies we plot the imaginary part of the susceptibility in Fig. [IJc). 
For Cr, Mn, Fe and Co the Heisenberg model gives -19.7 meV, 15.4 meV, 39.2 meV and 
-33.1 meV for the optical mode frequencies respectively. 

To conclude, we have shown that a simple approach, based on TD-DFT and the KKR-GF 
method, can be used to extract dynamics magnetic susceptibilities. This eliminates the need 
to acquire the parameters which enter empirical tight binding approaches, since a description 
of the electronic structure is embedded into our single particle Green functions. The size of 
matrices in the Dyson equation are small enough to permit calculations of large nanostruc- 
tures in future applications. We have developed an identity that leads to a numerically stable 

8 



method of extracting the proper Coulomb interaction to be used in the Dyson equation from 
which the full dynamical susceptibility is obtained. As an application, 3d adatoms and dimers 
deposited on Cu(OOl) surface were investigated from first-principles. 
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FIG. 1: In (a) is shown the variation of the product of U with the different values of the moments. 

n 

The values obtained for the adatoms fall on the diagonal as estimated in the curve of Himpsel|15l|. 
In (b) the imaginary part of X is plotted for every adatom when applying an additional magnetic 
field along the z-direction corresponding to a Larmor frequency of 13.6 meV (dashed line). In (c) 
are shown Imx calculated for the four dimers. The optical modes, estimated for Mn and Fe from 
a Heisenberg model, are represented as dashed lines. 
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